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Abstract 


Methods for predicting the onset of flutter during an experiment are traditionally applied treating 
the data as deterministic values. Uncertainty and variation in the data is often glossed over by using 
best-fit curves to represent the information. This paper applies stochastic treatments to wind tunnel 
data obtained for the Piezoelectric Aeroelastic Response Tailoring Investigation model. These 
methods include modal amplitude tracking, modal frequency tracking and several applications of the 
flutter margin method. The flutter margin method was developed by Zimmerman and 
Weissenburger, and extended by Poirel, Dunn and Porter to incorporate uncertainty. Much of the 
current work follows the future work recommendations of Poirel, Dunn and Porter. 
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Nomenclature 

modal amplitude 
flutter margin 
lower bound 

number of degrees of freedom (number of experiments) 
probability of event occurring 
sample standard deviation 
sample variance 

value of the random variable associated with a standardized normal distribution function 

standardized Gaussian random variable 

upper bound 

random variable 

mean value 

significance level 

modal decay rate, ~Ca) n 

number of standard deviations from the mean 

uncertainty on the mean for a Gaussian process 

standard deviation for a Gaussian process 

shape parameter for a gamma distribution 

scale parameter for a gamma distribution 

generic random variable 

modal frequency 

chi-squared distribution 

modal damping 


Subscripted variables: 

fx (**0 = probability density function of the random variable, X 

u a/2 = value of the random variable associated with a standaridized normal distribution function; 

the probability of a random variable, U, drawn from the standardized normal distribution 
assuming a value less than u a/2 is 100*(l-a/2)%. 


1 Research Engineer, Aeroelasticity Branch, Mail Stop 340, AIAA member. 
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Pb P2 

a UB 

COi, ©2 

2 

Zn-l,a/2 


decay rate of modes 1 and 2 

upper bound on the standard deviation confidence interval 
frequency of modes 1 and 2 

value of the chi- squared distribution where the probability of a random variable being 
greater than the distribution value is a/2, for n-1 degrees of freedom 


I. Introduction 


E XPERIMENTALLY determining the flutter speed of a flight vehicle or wind tunnel model can be 
a harrowing experience, as the onset of flutter can be rapid and destructive. To reduce the risks to the 
pilot, vehicle, model and facility, many methods of predicting the onset of flutter are used by test engineers. 
Traditionally, subcritical data and the flutter onset condition have been treated as deterministic values, 
although the experimental measurements may contain uncertainties due to data acquisition or interpretation 
and variability due to changes in the conditions or properties. 1,2 

Modal parameter uncertainties and their influence on flutter onset prediction were examined by Poirel, 
Dunn and Porter in 2005. 3 They extended the flutter margin method developed by Zimmerman and 
Weissenburger, 4 examining the effects of uncertainty in measured modal data. The current work parallels 
the application of the statistical methods demonstrated by Poirel, Dunn and Porter and attempts to address 
their recommended areas for future work. These areas include applying the methodology to a wind tunnel 
model with experimentally determined flutter speed, and exploring the concept of confidence intervals. 
Additionally, the overall framework they developed for incorporating uncertainties into flutter onset 
prediction is applied to inverse amplitude tracking. Theoretical stochastic models for the flutter dynamic 
pressure probability density function are examined, focusing on the Gaussian and gamma distributions. 
Statistical confidence bounds are computed for the Gaussian model and their implications are discussed. 

The current paper distinguishes itself from the prior work by applying the methods to experimental data, 
where flutter conditions have been measured. Additionally, confidence intervals are determined and their 
relevance to flutter predictions is assessed. 

II. The Piezoelectric Aeroelastic Response Tailoring Investigation 


This work utilizes data obtained in the Transonic Dynamics Tunnel (TDT) at NASA Langley Research 
Center in March and November 1994. 5 The Piezoelectric Aeroelastic Response Tailoring Investigation 
(PARTI) model was designed as a testbed to demonstrate flutter suppression using piezoelectric actuators. 
The model was constructed with a graphite-epoxy composite sandwich plate as the primary load-carrying 
structure. Six fiberglass aerodynamic shell elements were bolted to the plate at single contact points near 
the leading and trailing edges, and spanned with foam inserts covered with latex. The model is shown in 
Figure 1, mounted in the TDT test section. 



Figure 1. PARTI model mounted in wind tunnel 
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Data from testing the PARTI wing lends itself well to stochastic treatment because significant variations in 
flutter condition and modal properties occurred. The most significant illustration of model change is that 
the flutter dynamic pressure dropped by 10% between the March and November tests, Table 1. A change in 
the natural frequencies, air off, had also been observed relative to data taken prior to the March test. The 
repeatability of the aerodynamic shell assembly was questioned, and it was also suspected that the 
composite plate exhibited creep over the course of the program. These potential sources of variation and 
change were investigated. Removal and reassembly of the model did not significantly change the natural 
frequencies or the flutter condition. Non-destructive evaluation of the primary wing structure revealed no 
damage. It was decided to proceed with the test and treat the differences from the first test to be permanent 
changes to the model. 

Table 1. Flutter characteristics for the two wind tunnel tests 


Flutter Parameters for 
Stowed Configuration 

First Test 
(March) 

Second Test 
(November) 

Dynamic Pressure (psf) 

85 

76 

Frequency (Hz) 

6.5 

7 


Because the flutter speed had changed so significantly from the prior test, periodic flutter passes were made 
throughout the test to determine if there were changes due to variability in the tunnel conditions, or 
fluctuating model properties. Nine experimental flutter points, obtained during the second test, are detailed 
in Table 2- the flutter points occurred between 73 and 78 psf. All data was obtained using air as the test 
medium at atmospheric pressure. 


Table 2. Measured flutter conditions from November test 


Run # 

Flutter dynamic 
pressure (psf) 

Notes 

2 

75+ 

+ Tunnel condition was changing as flutter occurred; 
dynamic pressure was greater than 75 psf 

4 

75+ 

Tunnel condition changing between 75 and 76 psf 

20 

76 

Flutter point at stabilized tunnel condition 

24 

73 

Two passes, both flutter at 73 psf stabilized condition 

28 

74 

Two passes, both flutter at 74 psf stabilized condition 

35 

75+ 

Tunnel condition changing between 75 and 76 psf 

36 

78 

With trip strips 


Although the flutter condition did not repeat exactly during each of these runs, the data indicated that a 
permanent change in the model had occurred. This required acquisition of new model identification data, 
and redesign of the control laws. This effort produced a large data set, with the engineers particularly 
focused on observing variation in the modal and flutter characteristics. The remainder of this paper focuses 
on open loop data acquired during the November test, and treating the data from a statistical standpoint. 

III. Inverse Amplitude Tracking Method 


A. On-line prediction 

The procedure for conducting a flutter pass during the PARTI testing was to start at a subcritical dynamic 
pressure and acquire data. The data was analyzed in near real time to produce modal amplitudes, 
frequencies and dampings. On-line flutter onset predictions were updated and the dynamic pressure 
increased. This was repeated until flutter was encountered. 

One on-line method for predicting flutter onset is to track and curve fit the inverse of the modal amplitudes. 
In the PARTI test, data was acquired from accelerometers; a Fourier transform was applied to the data. 
Peaks appearing in the data correspond to the response at the aeroelastic frequencies of the system. The 
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peak value at each of the aeroelastic modal frequencies was inverted and plotted. Flutter can be indicated 
by extremely large modal amplitudes or, correspondingly, it can be indicated by near-zero inverse modal 
amplitudes. Figure 2 shows a data set obtained from a single flutter pass. 

The inverse of each amplitude for the first two aeroelastic modes is plotted versus dynamic pressure. 

Flutter onset predictions were generated by linear fits to each mode and extrapolated to the x-axis crossing. 
Throughout the test, it was found that the first mode’s amplitude was the most reliable for flutter 
prediction; this is illustrated by the on-line prediction results presented in Table 3. The predictions shown 
in the table were generated by linear fits using all of the data at and below the corresponding dynamic 
pressure. The experimentally determined flutter condition for this run was between 75 and 76 psf. Using 
all of the data up to 75 psf produces a slightly non-conservative flutter prediction, 76.2 psf. 

Examination of mode 2 illustrates an important point regarding the use of the inverse amplitude tracking 
method- engineers rarely rely on an automated linear fit to the full data set. The data at 30 psf skews the fit 
to data that is otherwise very linear. Typically, the engineer monitoring the process would ignore the 30 psf 
data point and provide a flutter estimate based on a modified linear fit. The revised prediction, however, is 
still non-conservative. 

The remainder of the discussion of the inverse amplitude tracking method will use only data corresponding 
to the first aeroelastic mode. 



Figure 2. On-line flutter prediction, inverse amplitude tracking, Run 2 


Table 3. On-line flutter predictions, Run 2 



Mode 1 

Mode 2 

Dynamic 

Pressure 

(psf) 

1 /amplitude 

Flutter 

prediction 

(psf) 

1 /amplitude 

Flutter 

prediction 

(psf) 

29.75 

100.0 


41.0 


51 

48.5 

71.0 

49.3 

-75.5 

60 

32.6 

73.5 

37.5 

1704.7 

71 

11.1 

75.1 

23.6 

150.6 

75 

5.6 

76.2 

18.9 

121.8 


B. Stochastic post-test data analysis: inverse amplitude histogram method 

While the inverse amplitude method described above worked well for the PARTI test, the data for each run 
had to be used in isolation. When multiple runs were compared, it was observed that there was significant 
scatter in the data. Normalized inverse amplitude data sets from eight flutter clearance passes are presented 
in Figure 3, illustrating this scatter. 
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Dynamic Pressure, psf 

Figure 3. Inverse amplitude tracking method for flutter clearance 

With all of the data from the test in hand, a stochastic treatment for the inverse amplitude data was 
developed and applied in a manner similar to the on-line prediction technique. Choosing one measured 
value of the amplitude at each of the test conditions, a linear fit was performed to determine a flutter onset 
prediction. This process was repeated, exhaustively choosing new values of amplitude at each of the test 
conditions and performing linear fits and predictions. The number of flutter dynamic pressure predictions 
is the number of degrees of freedom in this experiment. 

A normalized histogram of the flutter onset predictions is shown in Figure 4. The histogram approximates 
a probability density function, (PDF). Statistics of the process were calculated, yielding a mean of 83.9 psf 
and a standard deviation of 3.97 psf. The median value, influenced by the presence of a right-side tail is 
83.34 psf. The peak indicates the statistical mode value, 82 psf. The distribution is reasonably 
approximated as Gaussian; this idealized distribution is also shown in Figure 4. The cumulative 
distribution function, (CDF), Figure 5, is obtained by integrating the area under the PDF. The graph shows 
the probability that flutter has not occurred at a dynamic pressure lower than a given x-axis value. 

Stochastic treatment of the amplitude data produces a non-conservative prediction of the flutter dynamic 
pressure. There are two ways to observe this from the data. The experimental CDF, Figure 5, shows the 
probability that flutter will occur below a given dynamic pressure. The percent probability that flutter will 
not occur below a given dynamic pressure is determined by subtracting the vertical axis value from 100%. 
The value of 5% probability shown on the vertical axis corresponds to 78 psf on the horizontal axis. 
Therefore, there is a 95% probability that the model has not fluttered at 75 psf. 

The second method of observing the flutter prediction is to consider the Gaussian model of the process. 
Using the statistics to provide 95% probability that the system has not fluttered requires that the lower 
bound be calculated by subtracting twice the standard deviation from the mean. This gives a lower bound 
of 76 psf as the flutter dynamic pressure. It is not conservative to say that 95% of the time, the system will 
not flutter at conditions below this- recall from Table 2 that out of all nine points where flutter occurred, 
only one fluttered above this condition. 

Poirel, Dunn and Porter suggest using the mode as the central statistic to induce conservatism into the 
process. Applying twice the standard deviation to the mode gives a 95% confidence lower flutter bound of 
74 psf, which is more conservative. 


Heeg_SDM2007.doc 


5 


0,15 



Figure 4. Histogram of flutter predictions, using Figure 5. Cumulative distribution function, inverse 
inverse amplitude tracking amplitude method, linear fit 

Statistical treatment of the inverse amplitude data produced results which were not conservative. There are 
two primary causes: 1) A linear fit, constrained to use randomly combined data at each subcritical 
condition, did not generally give a conservative flutter prediction; engineering judgment generally needs to 
be applied to inverse amplitude data to accurately assess the flutter condition; 2) A Gaussian distribution 
can not fit the data accurately due to the skewness of the distribution of the flutter dynamic pressures. 

IV. Frequency Monitoring 

Modal frequency tracking was used as a qualitative guide for predicting the flutter condition during the 
PARTI test. Flutter onset was predicted by visual extrapolation of the frequencies to their coalescence. 

Knowledge of the analytical trends or the behavior on previous flutter passes was essential to accurately 
using this information. On-line data acquired for monitoring the frequencies of the first three aeroelastic 
modes is shown in Figure 6. This data corresponds to the amplitude data shown previously in Figure 3. 



Dynamic Pressure, psf 

Figure 6. Flutter clearance run, on-line data analysis tracking the frequencies of the first three 

modes 

This qualitative method can not be adapted to use stochastic techniques, but provides the data essential for 
using the flutter margin method. All of the natural frequency measurements obtained during the test for the 
first 2 modes are shown in Figure 7, illustrating a significant amount of scatter in the measurements. 
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Figure 7. Scatter in modal frequency measurements 


Histograms of the frequency information acquired at 50 psf are presented in Figure 8. Although neither 
distribution is precisely Gaussian, the statistics of the process were computed using this assumption. The 
mean, variance, mode and median values for each of the modal frequencies are given in Table 4. 




Frequency of Mode 1 (Hz) Frequency of Mode 2, Hz 


a) Mode 1 b) Mode 2 

Figure 8. Frequency histograms at 50 psf 


Table 4. Statistical Properties of Freq 

uencies at 50 psf 



Mean 

Variance 

Median 

CO n 

Mode 1 

3.868 

0.017 

3.875 

Mode 2 

11.107 

0.059 

11.0125 


The statistics were computed for the natural frequencies at each of the test conditions from 50 to 75 psf. 
Using the mean and the standard deviation, 95% probability values (+/- 2 a) on the frequencies were 
computed and are shown in Figure 9. This figure provides an interpretation of the potential variation in the 
flutter dynamic pressure. Using the 95% probability values, and assuming the worst-case flutter scenario, 
the dashed lines nearest each other represent the frequency characteristics of the aeroelastic system. Using 
these trends, flutter is imminent at 75 psf. However, if the least conservative case is assumed, the outer 
dashed lines represent the frequency characteristics of the aeroelastic system. While they have started to 
migrate towards each other, it is not evident that the onset of flutter will occur prior to 80 psf. 

The consequences of assuming that the Gaussian statistics accurately represent the experimental data can 
be seen in Figure 9b. Here, the data from Figure 7, which shows the experimental variation in the 
frequencies, has been overlaid with the Gaussian 95% probability bounds. Where the experimental data 
points lie outside the bounds, it is obvious that the Gaussian representation is not accurate- representing the 
frequency data by the statistics of a Gaussian process is not a good methodology to employ. 
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Figure 9. Statistical calculations of frequency versus test condition 


V. Flutter Margin Method 


Zimmerman and Weissenburger 4 defined a flutter stability parameter which they demonstrated to be highly 
effective at predicting flutter onset for classical bending-torsion flutter. Their method was developed based 
on an idealized bending-torsion system subject to quasi-steady aerodynamics. They developed a parameter 
that they termed the flutter margin, which incorporated measured frequencies and decay rates at a given 
dynamic pressure. Using the flutter margin at a series of dynamic pressures, they extrapolated to find the 
value of dynamic pressure that had a flutter margin of 0, and thus produced a prediction of the flutter 
dynamic pressure. They also showed that the flutter margin was a second order function of dynamic 
pressure. At a given test condition, the flutter margin is calculated using Equation 1. 


F = 


2 2 
CD 2 ~ CD\ 


{ Pi -Pi 


+ 4 PxPx 


' 2 2 A 

CO 2 + (D\ 


+ 2 


Pi + P\ 


Pi P\ 
Pi + P\ 


V 2 1 N f 

0 ) 2 - CO\ 


\ f \ 

P 2 +Pl 


Equation 1 


Poirel, Dunn and Porter recognized that the flutter margin method was sensitive to scatter in the measured 
values of the frequencies and the decay rates. They applied the flutter margin method in a manner that 
captured uncertainties in the measurements. 

They utilized two basic methods to generate flutter predictions: 1) A of values for the modal parameters at 
each test condition was used to compute the flutter margin; the mean value of flutter margin at each 
dynamic pressure was used to extrapolate to the flutter prediction; 2) The range of values for the modal 
parameters at each test condition was used to compute all possible values of flutter margin; all possible 
combinations of flutter margin at each of the test conditions were then used to produce an extrapolated 
flutter prediction. An identical application of these methods follows, using measurements from the PARTI 
experiment. 
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A. Mean flutter margin method 

The first method requires determining the mean value of the flutter margin at each dynamic pressure. 

These mean values are then plotted as a function of dynamic pressure and a quadratic curve fit is applied. 
The extrapolation of this curve fit to 0 produces the estimate of the flutter dynamic pressure. 

As previously shown, numerous measurements of frequency were obtained at each of the test conditions. 
Modal damping was measured only sporadically during the PARTI test. Fortunately, Poirel, Dunn and 
Porter determined that the damping variation has much less impact on the flutter margin calculations than 
the frequency variation. Thus, the mean values of the measured modal damping at each dynamic pressure 
were used in computing flutter margins. 

For each test condition, all of the combinations of measured frequencies at that condition were used to 
compute values of flutter margin. This gave rise to a large number of possible flutter margin values. The 
mean of these values at each test condition was then curve fit and a flutter prediction was extrapolated. The 
results are shown in Figure 10. The fit produces a very good estimate of the flutter condition. This is true 
even when only the lowest 2 dynamic pressures are used in the curve fit. 



Figure 10. Mean flutter margin extrapolation results 


B. Flutter margin histogram method 

The flutter margin histogram procedure is similar to the inverse amplitude histogram procedure previously 
described. Applying the second method begins the same as the first method. All possible values of flutter 
margin at each condition were generated, exhaustively using the measured frequencies. The damping at 
each test condition was taken to be the mean measured value. Example histograms of the flutter margin 
results are shown in Figure 1 1 . Figure 11a shows all calculated flutter margin values at 50 psf; Figure lib 
shows all values for 70 psf. 



Figure 11. Histogram of flutter margin at example test conditions 
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At this point, the second method deviates from the first. At each of the test conditions, there are large 
numbers of values for the flutter margin. In this analysis, one value of the flutter margin was selected at 
each of the test conditions. The value of flutter margin at each test condition was selected by random draw 
from all of the calculated values. A quadratic polynomial was fit to these flutter margin values, as a 
function of dynamic pressure to produce a prediction of flutter dynamic pressure. This process was 
repeated to produce 100,243 possible flutter values. The number of flutter dynamic pressure values is the 
number of degrees of freedom in this experiment. A normalized histogram of the results is shown in Figure 
12 . 

Identical processes were applied to the data, using only information available up to 65 psf and then up to 70 
psf. The statistics of these processes are given in Table 5. Of the three central statistics, the mode is 
always the lowest value. The mean and variance values vary substantially as more test conditions are 
included in the analysis. As data is included closer to the flutter condition, the variance in the predictions 
decreases. This is expected, because the subcritical data at 75 psf indicates that flutter will occur 
imminently. However, including subcritical data at 75 psf all but precludes consideration that flutter could 
already have occurred. Recall from the flutter points discussed earlier that numerous flutter points did 
occur below 75 psf. 



Dynamic pressure 

Figure 12. Histogram of flutter dynamic pressure predictions, Flutter margin histogram method 


Table 5. Statistics of predicted flutter dynamic pressures (psf) 


Data Included 

Mean 

Variance 

Median 

Mode 

50-65 psf 

114.3 

1.2e7 

69.9 

66.5 

50-70 psf 

76.9 

17.5 

75.8 

74.2 

50-75 psf 

77.1 

3.4 

76.6 

76.0 


VI. Fitting a Theoretical Model to the Flutter Dynamic Pressure 

Calculating the mean and variance assumes an underlying theoretical Gaussian model of the distribution. 
A Gaussian process has symmetric distribution about the mean, with the sharpness of the peak defined by 
the variance. Using the estimates for the mean and the variance, the theoretical Gaussian distribution is 
plotted with the flutter dynamic pressure histogram in Figure 12. It shows the shortcomings of this model; 
the Gaussian process captures neither the peak nor the skewness of the distribution. 

Non-symmetric distributions have been used for modeling this type of process. A gamma distribution 
resembles the shape of the flutter dynamic pressure histogram. However, obtaining a reasonable fit to the 
data required that a threshold value or a bias be applied in order to avoid numerical issues associated with 
calculating the gamma function for large values of shape parameter. A maximum likelihood, (MLE), 
parameter estimation method 6 for the gamma distribution parameters 7 was programmed. Using the biased 
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data, maximum likelihood estimates of the shape and scale parameters were 6.971 and 1.673, respectively. 


The plot in Figure 13 shows the histogram of the biased data, using the solid blue line with circles. The 
MLE gamma function is shown by the dashed cyan line. Note that the peak value is approximately 75% of 
the peak value from the histogram of flutter dynamic pressure. It is also offset from the peak, and is 
broader than the histogram. A second gamma distribution, which was manually generated and is not 
optimized in any sense other than appearance, is shown by the solid magenta line. The shape parameter is 
9.0; the scale parameter is 2.7. The peak value is much closer to that of the histogram, but note how the 
distribution has a left shift. This curve will give a more conservative estimate of the flutter dynamic 
pressure than the MLE result. 


Shape 

Parameter 

n 

Scale 

Parameter 

X 

Generated 

through 

Gamma Distribution 

6.971 

1.673 

Maximum Likelihood 

Gamma Distribution 

9.0 

2,7 

Manual Fit 

■©■ Histogram of Data 






73 77 81 85 89 93 


Flutter dynamic pressure, psf 

Figure 13. Theoretical model fit to flutter dynamic pressure histogram 


A further test of the goodness of a theoretical model is provided by examining a quantile- quantile plot. 
From the NIST/SEMATECH e-Handbook of Statistical Methods 8 , “The quantile-quantile (q-q) plot is a 
graphical technique for determining if two data sets come from populations with a common distribution. A 
q-q plot is a plot of the quantiles of the first data set against the quantiles of the second data set. By a 
quantile, we mean the fraction (or percent) of points below the given value. That is, the 0.3 (or 30%) 
quantile is the point at which 30% percent of the data fall below and 70% fall above that value.” 

A q-q plot is presented in Figure 14, where there experimental data is compared to itself, the experimental 
data is compared to the maximum likelihood estimate of the gamma distribution, and the experimental data 
is compared to the manual fit gamma distribution. The experimental data compared to itself shows what 
the data would look like if one of the theoretical models were an exact fit to the data. 

When a quantile plot goes above the exact fit line, the theoretical model is non-conservative. The black 
dashed lines in Figure 14 give examples for the two theoretical distributions. From the vertical axis, the 
horizontal line at 0.6 intersects the maximum likelihood distribution at a point corresponding to 0.5. This 
says that 60% of the data points, or values for flutter dynamic pressure, are below the dynamic pressure 
where the model estimates that only 50% of them he below. Thus, the maximum likelihood model is not 
conservative, even though the overall fit to the data appears better than the manual fit. 

When the quantile plot goes below the exact fit line, the theoretical model is conservative, as in the case of 
the manual fit gamma distribution. From the vertical axis, the horizontal line at 0.35 intersects the manual 
fit distribution at a point corresponding to 0.5. This says that 35% of the data points, or values for flutter 
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dynamic pressure, are below the dynamic pressure where the model estimates that 50% of the values would 
lie below. Thus, the manual fit model gives a conservative estimate of the flutter dynamic pressure 
distribution. 
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X 
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Figure 14. Quantile plot using gamma distributions 

VII. Probability levels associated with flutter 

The quadratic curve fits to the flutter margin data produce some cases which do not flutter. Table 6 shows 
the results for the three ranges of test conditions which were considered. Using only data up to 65 psf, only 
38.5% of the curve fits indicated that the system would flutter. As more data was included, the probability 
of flutter occurring rose to 99.1%. While this indicates that including higher dynamic pressure information 
in the subcritical data set produces a better prediction of flutter, closer scrutiny is required. 

The data was also examined to determine the dynamic pressure bound corresponding to 95% probability of 
being flutter-free. This was done for each of the three ranges of test condition, corresponding to the three 
data rows in Table 6. The 3 rd column of Table 6 shows the results obtained when the information in the 
experimental CDF is used. As more data is incorporated, the dynamic pressure bound increases. It is 
interesting to note that each row of data indicates that the flutter dynamic pressure is very close to the 
highest condition considered. 

The final two columns of the table give the dynamic pressures that correspond to 95% probability of being 
flutter-free if the statistics of the estimated Gaussian processes are used. The high variance using the first 
range of test conditions produces a negative dynamic pressure lower bound. The second data set considers 
conditions in the subcritical data set that lie strictly below the lowest obtained experimental dynamic 
pressure. The results produce conservative estimates for the bound on flutter: 68.5 psf using the mean 
value and 65.9 psf using the mode value. The third data set, considering all subcritical data available, gives 
higher values of the 95% probability lower bound. 


Table 6. Probabilities associated with flutter conditions 




Dynamic pressure corresponding to 
95% probability of being flutter-free, (psf) 

Test conditions 
included 

Probability that 
flutter will occur, 
(% of cases that 
fluttered) 

Using CDF 
data directly 

Based on 
Gaussian 
process, 
Mean - 2a 

Mode - 2a 

50-65 psf 

38.5 

66.5 

N/A 

N/A 

50-70 psf 

90.2 

72.4 

68.5 

65.9 

50-75 psf 

99.1 

75.3 

76.0 

72.3 
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Another way to examine the probability, is to determine dynamic pressures corresponding to different 
probability levels for flutter occurrence. This has been done using the gamma distributions, and using all 
subcritical data available. The CDFs of the data and the gamma distributions are shown in Figure 15. 
Example calculations for the conditions corresponding to 95% probability of being flutter- free are 
illustrated in the close-up, part b of this figure. 


The dashed line originating on the vertical axis shows when the CDFs indicate that there is only a 5% 
probability of flutter occurring. This corresponds to a 95% probability of being flutter- free. The first CDF 
that this line intersects is the manual fit gamma distribution, 74.7 psf. Next, it intersects the maximum 
likelihood gamma distribution at 74.9 psf. Finally, it intersects the CDF produced directly by the quadratic 
curve fits to the data, 75.3 psf. This progression indicates that both gamma distribution models produce 
conservative values for the 95% flutter- free estimates. 
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a) Showing full range of possible flutter conditions b) 95% Probabability conditions 

Figure 15. Cumulative distribution functions for the gamma distributions 


However, reviewing the quantile plot, it is apparent that this will not be true for the maximum likelihood 
model at other probability levels. In many statistical analyses, it is customary to show data corresponding 
to different a-levels. Although these are strictly applicable only to a Gaussian process, these levels 
correspond to levels of probability. Using 1-, 2- and 3 -a levels, Table 7 summarizes the dynamic 
pressures corresponding to the associated levels of probability. The columns show the results for direct use 
of the flutter margin data CDF, the Gaussian model and the two gamma distributions. 


Notice that for each probability level, even though the Gaussian model is the poorest representation of the 
data, the symmetric properties of the distribution allow it to produce the most conservative values for flutter 
onset. 


Table 7. Flutter predictions corresponding to probability levels, using data up to 75 psf 




Dynamic Pressure (psf) 




Theoretical Model 

Probability of 
being flutter- 
free 
(%) 

Sigma level 
of a Gaussian 
distribution, 
(for reference 
only) 

Flutter margin 
histogram 
results, directly 

Gaussian, 
point estimates 

Gamma, 
maximum 
likelihood fit 

Gamma, 
manual fit 

38 

72 

77.02 

76.21 

77.42 

76.5 

68 

1 

76.11 

75.29 

76.28 

75.70 

95 

2 

75.27 

73.46 

74.92 

74.70 

99.7 

3 

74.68 

71.63 

74.06 

74.04 
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Also note that the maximum likelihood method for determining the shape and scale parameters did not 
build conservatism into the model, except at the very high and very low probability of flutter occurring. 
Manual adjustment of the parameters allowed this conservatism to happen. In a case like this, where 
conservatism is essential, a different fit criterion might be a better choice. The gamma distribution may 
provide the best fit that can be achieved with a 3 -parameter model, where the three parameters are the 
shape, scale and threshold. Higher order distributions, as suggested by Pimentel et al, 9 Bandourian et al, 10 
and Clark et al. 11 would likely provide better fits to the data. 

VIII. Confidence Intervals 


The threshold gamma distribution model fits the flutter margin data much better than the Gaussian 
distribution model. However, the equations for calculating the confidence intervals for a gamma 
distribution are not widely-published. For the first cut at determining confidence levels for the flutter 
dynamic pressure, it was assumed that the distribution could be estimated as a Gaussian process. Methods 
for determining confidence intervals for Gaussian processes are well known and easily accessible. 6,8,12,13 

There are several approaches to estimating statistical parameters from a set of sample measurements. The 
first approach generates a single value to estimate the parameter. Here, the sample mean has been 
calculated as an estimate of the process mean and the sample variance has been calculated as an estimate of 
the process variance. 

A second approach to estimating a parameter is interval estimation, which produces bounds on the 
parameter value. These bounds give numerical information for the parameter and also indicate the level of 
confidence that can be placed on the numerical value based on the samples. The confidence interval on the 
mean value and on the variance can be computed using the normal distribution function and the chi-squared 
distribution function, respectively. Upper and lower limits of an interval are defined such that, with a 
specified level of confidence, the true parameter will lie in this interval. The interval (LB,UB) is the 100(1- 
a)% confidence interval for a parameter 0, 6 as expressed in Equation 2. 

P(LB < 0 < UB) = 1- a Equation 2 


In evaluating the flutter points, the goal is to determine the lowest value at which flutter is likely to occur 
with a given probability. Thus, a lower bound on the flutter dynamic pressure is calculated using the 
general form of Equation 3. 


LB = X -K-y* (Jjjb 


Equation 3 


A graphical illustration applying these confidence intervals is illustrated in Figure 16. The plot shows a 
standardized Gaussian distribution probability density function, f x (x) with a solid blue line. The solid red 
line at X=0 shows the mean value, X , of the process. The dashed red line, just to its left indicates the 
uncertainty associated with the mean of the process, k. It is located at X - k . The solid green lines 
indicate the X -la values, which correspond to the 68% confidence interval if the statistics of the process 
are known. The dashed green line to the left of the solid green line is a lower bound for 68% confidence 
which includes the confidence interval on the mean. The dashed cyan line further to the left then shows the 
influence of applying the confidence interval of the variance instead of assuming that the computed 
variance is known with 100% confidence. The upper bound on the variance is used because that makes the 
variance a larger number, representing less certainty. 
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Figure 16. Confidence intervals on a normal distribution 

The values for k and a are determined in the following two sections. 

A. Confidence interval on the mean 

Estimating the confidence interval for the mean of a normal distribution with a known variance uses the 
symmetry properties of the normal distribution and assumes use of a minimum width interval. The 
uncertainty on the mean is computed by transforming the process to a standardized random variable, 
assuming the maximum likelihood estimator for the mean is itself normal, and the probability density 
function of a Gaussian process 6 . This results in the uncertainty on the mean as expressed in Equation 4. 

The confidence interval bounds are determined by the uncertainty subtracted and added to the mean. 

_u a i 1 G Equation 4 

Note that this expression is strictly valid only for processes where the variance is known. The student’s t- 
distribution and the sample variance are employed otherwise. Because there are a large number of degrees 
of freedom in this problem, the numerical results were identical to the results using the known variance. 6 

B. Confidence interval on the variance 

The confidence interval for the variance is determined in a similar manner, but utilizes the chi-squared 
distribution function. Upper and lower bounds define the 100(l-a)% confidence interval for the variance, 6 
Equation 5. 


(*-l)S 2 

Zn-l,a/2 


2 («-l)S 2 

<<J < 


1-a 


Zn-l,l-a/2 


Equation 5 


The upper bound on the confidence interval for the variance is taken from Equation 5 and shown in 
Equation 6. 
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Equation 6 


M 2 


(n-l)s 2 

2 

Zn-l,l-a / 2 


This study contained over 100,000 degrees of freedom in predicting the flutter dynamic pressure. As the 
number of degrees of freedom increases, the bounds squeeze in towards the point estimation of the 
variance. For an infinite number of degrees of freedom, the bounds converge to the point estimate. This 
corresponds to the situation where the variance of the process is known, and no uncertainty interval on the 
variance is required. Also, the chi- squared distribution becomes a Gaussian distribution as the number of 
samples goes to infinity, based on the central limit theorem. 6 Equation 7 shows the behavior of the upper 
bound on the confidence interval as the number of degrees of freedom becomes large. 


, ( *~ l) -> l, so 

% n-l,l-a / 2 


Equation 7 


C. Gaussian confidence results 

Table 8 shows results utilizing point estimates of the mean and variance and incorporating the confidence 
intervals into the data. Each row in the table shows the lower bound on the flutter dynamic pressure for the 
percent confidence given in the second column. The fifth column shows the results when point estimates 
of the mean and variance are used; the sixth column shows the results when the confidence intervals on the 
mean and variance are included in the computation. For illustration, 100 degrees of freedom were assumed 
in calculating the bounds on the variance. Note that incorporating these additional uncertainties gives 
slightly more conservative flutter onset predictions. 

Table 8. Uncertainty parameters and lower bound on flutter dynamic pressure; Gaussian 
theoretical model 


Sigma 

Level 

(Y) 

% 

Confidence 

Uncertainty 
on the mean 
(K) 

& UpperBound 

Lower Bound on flutter dynamic pressure 
(psf) 

Using point 
estimates of mean 
and variance 

Incorporating 
confidence 
intervals on the 
mean and variance 

Vi 

38 

0.0029 

3.62 

76.21 

76.17 

1 

68 

0.0058 

3.9 

75.29 

75.14 

2 

95 

0.0113 

4.51 

73.46 

72.86 

3 

99.7 

0.0168 

5.41 

71.63 

70.13 


IX. Discussion 

The flutter margin method provides a quantitative method for utilizing frequency coalescence information 
available at subcritical test conditions. This work demonstrated that it is accurate for predicting the onset 
of classical coalescence flutter. This is true for application of the method both in a deterministic and a 
stochastic manner. 

A fundamental question to consider is whether or not stochastic methods can be used to conduct flutter 
testing. 

The stochastic methods provided good estimates of the flutter dynamic pressure, although slightly non- 
conservative. The non-conservatism is partially due to non- conservative modal trends, rather than the 
stochastic methods themselves. Applying these methods in an automated fashion removed the influence of 
engineering insight from data selection and interpretation. In real-time flutter testing, it is often this insight 
that builds the conservatism into the predictions. 
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Traditionally, when conducting flutter testing in a wind tunnel, you do not have a large amount of data 
available at the subcritical conditions prior to proceeding to flutter. Using stochastic methods would 
require a change in test procedure, acquiring numerous sets of modal characteristics at each subcritical 
condition. Another possible source of multiple data sets is to utilize multiple sensors at each subcritical 
condition. This eases the requirement to repeat the same test condition, but the near real time processing 
requirements for each test condition increase. 

For flutter testing where a Mach number-dynamic pressure boundary of flutter conditions are to be 
determined, stochastic methods may be suitable for determining the significance of the Mach number. If 
the subcritical data indicates that Mach number isn’t a significant factor, the flutter dynamic pressure would 
be invariant. This statistical methodology could thus indicate the onset of compressible effects. 

Using probability distribution models is a convenient approach because the stochastic behavior can be 
accomplished with only 2 or 3 parameters. However, the fits to the data are often not good. In this case, 
maximum likelihood parameter fits produced non-conservative gamma distribution models of the flutter 
dynamic pressure. To build conservatism into the model, a manual fit or another method of determining 
the coefficients needs to be used. 

Calculating the Gaussian confidence intervals demonstrated that for large data sets, the statistical 
confidence in the distribution parameters is high, even when the model does not accurately represent the 
data. Meeker and Escobar 2 state that uncertainty results from limited knowledge or limited data set. They 
separate it into uncertainty in model form and uncertainty in model parameters. “Uncertainty leads to 
uncertainty bounds- statistical confidence bounds quantify parameter uncertainty, but not model 
uncertainty.” In this case, the confidence intervals didn’t lend much insight to the problem. Rather than 
looking to the confidence intervals to determine the likelihood of flutter occurring, a measure of the 
goodness of the fit should be used to provide this insight. 

A noticeable shortcoming of the treatments presented here is that discrete variables were not incorporated. 
This is important in considering the range of conditions to incorporate into the subcritical data set. While 
incorporating all available subcritical data seemed like a good idea at the time, it unduly biased the results, 
forcing the dynamic pressure estimates to be high. By definition, subcritical data occurs below the flutter 
condition. By including the 75 psf data in all of the processes, there were only a minimal number of cases 
which produced flutter conditions below this value. A better approach would have incorporated knowledge 
of the number of times that the 75 psf data was post- flutter. 

X. Conclusions 

Experimental prediction of flutter onset using stochastic methods require a change in the process employed 
to conduct flutter testing. The methods evaluated here rely on numerous data points to quantify the 
uncertainty in modal characteristics at subcritical conditions. These methods are easily applied to in-hand 
historical data sets. 

Stochastic treatment of all of the data from the PARTI test illustrated the extent to which engineering 
judgment was relied upon to build conservatism into the flutter predictions. While sample automated fits 
shown in the paper provided good predictions of the flutter condition, these predictions were not 
conservative. 

One source of the non-conservatism is the definition and treatment of what constitutes the subcritical data 
set. Incorporating all available subcritical data forced the dynamic pressure estimates higher. A better 
approach would have incorporated knowledge of the number of times that the 75 psf data was post- flutter. 
Use of a discrete variable, indicating flutter or no flutter, might be useful for accomplishing this. 

The theoretical distribution of the flutter dynamic pressure requires additional consideration. The goodness 
of the theoretical distribution is a major factor in whether the overall prediction of flutter onset is 
conservative. While the threshold gamma distribution parameters can be manually tuned to provide a 
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conservative model, it is desirable to develop a parameter estimation process which includes this 
requirement. 

Using stochastic methods, it was demonstrated that probability levels can be associated with different 
flutter speeds, both neglecting and accounting for the confidence intervals. 

The methodology and conclusions of Poirel, Dunn and Porter are supported by analysis of the PARTI data. 
The methodology demonstrated an approach to handling uncertainty present in the data, in a manner that 
produced rational trends in the prediction of flutter onset. Both flutter margin methods were shown to 
produce good predictions of the flutter condition. The sensitivity of the quadratic fit extrapolation was also 
demonstrated by the percentage of cases that did not flutter. While the suggestion for development of 
confidence intervals was followed, the relevance of the confidence intervals for large data sets was shown 
to be negligibly small. 
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